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Abstract. We explore the possibility of putting constraints on dark energy models 
with statistical property of large scale structure in the non-linear region. In 
particular, we investigate the w dependence of non-Gaussianity of the smoothed density 
distribution generated by the nonlinear dynamics. In order to follow the non-linear 
evolution of the density fluctuations, we apply N-body simulations based on P^M 
scheme. We show that the relative difference between non-Gaussianity of w = —0.8 
model and that oi w = —1.0 model is 0.67% (skewness) and 1.2% (kurtosis) for 
R = 8h~^ Mpc. We also calculate the correspondent quantities for R = 2h~^ Mpc, 
3.0% (skewness) and 4.5% (kurtosis), and the difference turn out to be greater, even 
though non-linearity in this scale is so strong that the complex physical processes 
about galaxy formation affect the galaxy distribution. From this, we can expect that 
the difference can be tested by all sky galaxy surveys with the help of mock catalogs 
created by selection functions, which suggests that non-Gaussianity of the density 
distribution potentially plays an important role for extracting information on dark 
energy. 
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1. Introduction 

The formation of large scale structure in the Universe is a central element of astrophysics. 
The key idea in understanding the formation of the structure in the Universe is 
gravitational instability [H El El H]. Roughly stated, if the material in the Universe 
distributed irregularly, then the overdense regions provide extra gravitational attraction, 
draw material toward them and become more overdense. At present, this simple picture 
dramatically succeeds in explaining what is observed and current attention is focused 
entirely on the details of the gravitational instability process which strongly depends on 
the cosmological models. 

One of the useful means to extract the detail of structure formation is the statistical 
property of observable quantities. Among various statistical quantities, the density 
probability distribution function (PDF) of the density contrast plays a very important 
role because it includes information of non-Gaussianity. In the standard picture, the 
primordial density fluctuations originated from quantum fluctuations are stretched to 
large co moving scales during the inflation phase and assumed to be random Gaussian [5] . 
It is, however, well known that even though the density fluctuations remain Gaussian 
in the linear regime, it significantly deviates from Gaussian once the non-linear stage is 
reached [HllTj. Since structure formation is competition between the cosmic expansion 
and gravity, the growth of the density contrast as well as the evolution of its statistical 
quantities are intimately tied to the cosmic expansion law of the Universe. 

On the other hand, according to recent observations for type la Supernovae [8], 
the expansion of the Universe is accelerating. Combining measurements of Cosmic 
Microwave Background Radiation [9] and recent galaxy redshift survey [lOj, we are 
forced to recognize the existence of a cosmological constant, or kind of dark energy 
whose value is almost the same order of magnitude as the present density of the 
Universe [11]. (For reviews about the cosmological constant problem, see [I2]-) From 
the phenomenological viewpoint, in order to clarify the property of dark energy, it is 
important to constrain the effective equation of state of dark energy w by observations. 
Even though we have obtained the constraint as w < —0.90 (95 % confidence limit 
assuming w > —1) [11] by combining WMAP data with other astronomical data, in 
order to pin down the value of ly, it is necessary to propose other independent and 
complementary methods. 

For this purpose, as stated at the beginning, we suppose the statistical property 
of large scale structure of the Universe is helpful. Actually, we were the first to 
analyze the w— dependence of non-Gaussianity of the PDF of the smoothed density 
contrast in [13j where we considered simple constant-w dark energy models with 
w = —0.5, —0.8, —1.0, —1.2. (The case for the cosmological constant {w = —1) had 
already been considered in such as [TH UHl HE].) There, in order to see the effects 
of the non-linear dynamics analytically, as a first step, we adopted Lagrangian linear 
perturbation [171 [HI [IS] which describes the evolution of the density fluctuation in the 
quasi- nonlinear region better than second order Eulerian perturbation [201 [21]. 
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As a natural extension to the previous one, in order to go into the strongly nonlinear 
region where the difference of non-Gaussianity depending on w is expected to be large, 
numerical simulations seem to be necessary. Even though the simulations including 
dark energy are important, because of computational expense of the simulations, the 
first simulations including dark energy other than the cosmological constant appeared 
only in 2003. By now the aim of N-body simulations in the presence of dark energy 
are limited to the mass function of dark matter halo [221 l23l |2lj as well as the matter 
power spectrum [25j. Therefore, in this paper, we analyze the w-dependence of non- 
Gaussianity of the PDF of the smoothed density contrast with N-body simulations to 
cover the strongly nonlinear region where the results of [13] remain invalid. 

This paper is organized as follows. In section [2l we briefly summarize the effects 
of dark energy on the structure formation scenario and define the quantities which are 
necessary to estimate non-Gaussianity of the density distribution in section[3l The set up 
and results of N-body simulations are described in Section [Hand Section [5], respectively. 
Section [6] is devoted to conclusions. 



2. Effects of dark energy 

Here, we briefly summarize the effects of dark energy. We consider almost the same 
situation as our previous paper [13] where we assumed the following facts; 

(a) As for matter components, we consider (cold dark) matter and dark energy, 
i.e. ptot = Pm + Pde, where ptot, Pm and pde are total energy density of matter, each 
component of energy density of matter and dark energy, respectively. 

(b) Dark energy interacts with matter only gravitationally, i.e. the energy 
conservation of matter and dark energy holds independently. 

(c) The effective equation of state of dark energy w is constant for simplicity, even 
though in principle, we can relax this condition. 

(d) Curvature of the Universe is negligible, i.e. (/C = 0). 
Then, the background Friedmann equation can be written as 

= Hi [^]^o«-' + fiDi?o«-'^"+'T , (1) 

where a is a scale factor of the Universe, if = d/a is its Hubble parameter. We have 
also defined the density parameters for matter and dark energy: 

"(m, DE) = -^JpP{m, DE) ■ (2) 

The subscripts denote that the corresponding quantities are evaluated at the present 
time. From Eq. ([1]), the effects of dark energy affect the background cosmic expansion, 
which depends on the effective equation of state of dark energy w. 

Since structure formation is described as competition between the cosmic expansion 
and gravity, the change of the cosmic expansion modifies the growth of the fluctuations. 
In the following, we summarize the results obtained in structure formation, limiting to 
the linear perturbation theory, for simplicity. 
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In the comoving coordinates, the Poisson equation is given by 

V^'l' = ^QmH'a^6m, (3) 

where $ is the gravitational potential in the comoving coordinates, 6m = {Pm — Pmb) I Pmb 
is the density fluctuation of matter, p^h is the background density of matter. 

Combined with the continuity equation and Euler equation, the growth equation 
of the matter perturbation can be obtained 

'L + 2H6m - ^H^^m^m = 0, (4) 

where the dot denotes a time derivative. 

Therefore the growth of large scale structure depends on the background cosmic 
expansion law and it provides one of the constraints to the property of dark energy. 

Strictly speaking, in dark energy models with a dynamical field, one might wonder 
whether perturbations of the field might affect the growth of structure. However, since 
calculations of the galaxy power spectrum in such models suggest that the field is 
essentially smooth on these scales [261 1211, we assume that the perturbations of the 
field do not impact structure on these scales. 



3. Non-Gaussianity of the density distribution 

In what follows, since we are usually interested in the density fluctuations of matter, 
we will omit the subscript m unless otherwise stated. In order to analyze the statistical 
properties of the density contrast, we introduce a one-point probability distribution 
function of the density fluctuation field P{5) (PDF of the density fluctuation) which 
denotes the probability of obtaining the value 5. If 5 is a random Gaussian field, the 
PDF of the density fluctuation is determined completely as 

where cr^ = ((5 — (5))^) is the dispersion and ( ) denotes the spacial average. 

At initial, 5 is often treated as a random Gaussian field as a result of the generic 
prediction of inflationary scenario [28]. In the linear perturbation theory, Gaussianity 
of 5 is completely conserved if we start with a Gaussian initial condition, because each 
Fourier mode evolves independently according to the growth law Eq. (jl]). 

Once nonlinear terms are considered, however, the PDF deviates from the initial 
Gaussian shape. The point is that b is constrained to have a value 5 > — 1, otherwise 
the energy density p would be negative [1]. The Gaussian distribution ([5]) always assigns 
no n- zero probability to the unphysical regions with 5 < — 1. The error is negligible when 
a is small because the probability of obtaining the value 5 < — 1 is then very small, but, 
as fluctuations enter the strongly nonlinear regime with cr ~ 1, this error must become 
so important that the Gaussian distribution is only a poor approximation to the real 
distribution. What happens is that, as fluctuations evolve, mode-coupling effects cause 
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the initial distribution to skew, generating a long tail at high 6 while they are bounded 
at 5 = — 1. 

If the PDF deviates from Gaussian distribution, the cumulants of the one-point 
PDF of the density fluctuation field whose orders are higher than two become nonzero. 
Especially, the third and fourth order cumulants, which are defined as < 6^ >c=< >, 
< 5^ >c=< S'^ > —Sa^ mean the display asymmetry and non-Gaussian degree of 
"peakiness" , respectively, for a given dispersion. 

Since it is known that the scaling < 5" >cOC o"^""^ holds for weakly non-linear 
region, if we consider the gravitational clustering from Gaussian initial conditions |30j . 
we introduce the following higher-order statistical quantities [HE]: 



skewness : 7 



kurtosis : n = — — 

' t6 



at 

The merit of adopting these definitions is, as stated above, that they are constants in 
extremely weakly nonlinear stage which are given by Eulerian linear and second-order 
perturbation theory [Tl[30]. For example, in Einstein-de Sitter Universe, the skewness 
and kurtosis are given by 

60712 2x 
' 1323 ^ ' 

Even though it is shown that the skewness in weakly nonlinear region does not 
deviate much from 34/7 in several dark energy models [3T1 [32l [33] . the region we can 
observe is nearer and the nonlinearity becomes stronger. Therefore, in the following, we 
concentrate on these quantities by analyzing the nonlinear evolution of the fluctuations. 

It is worth noting that in the previous paper [13], we used other definitions for 
the skewness and kurtosis which are based on [6] in which < 5" >c is divided by a". 
The ones in the previous definition should be divided by a and cr^ for the skewness and 
kurtosis, respectively, to obtain the ones in the present definition. 



4. N-body simulations 

In our previous paper [1^, in order to analyze the non-linear dynamics, we applied 
Lagrangian perturbation theory where we can extract the quasi-nonlinear nature even 
if we consider only the linear order. Since we could obtain the naive w-dependence of 
non-Gaussianity of the density distribution in the weakly nonlinear region, as a natural 
extension, in this paper, we apply cosmological N-body simulations in order to obtain 
quantitatively plausible predictions. 

Below we describe the N-body codes, initial conditions, parameter choices and 
technical details of simulations. The results themselves are reported in Section [5l 

The numerical algorithm is applied by particle-particle particle-mesh (P^M) 
method whose code had been originally written by Bertschinger [3l]. The initial 
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Table 1. The relation between w and a%. This relation is derived from Eq. (fTTj). 



w 




-0.5 


0.56967 


-0.8 


0.663688 


-0.9 


0.698358 


-1.2 


0.813615 



Gaussian fluctuation is generated by COSMICS [35]. Until the quasi-nonlinear regime, 
COSMICS evolves the fluctuation within the scheme of Lagrangian linear perturbation. 
Then, we start N-body simulations at z ~ 30. 

For N-body simulations, we set the parameters as follows: 

Number of particles : = 128^ , 

Box size : L = 128/i"^Mpc (at a = 1) , 

Softening length : e = bOh^^pc (at a = 1) . 

The cosmological parameters at the present time (a = 1) are given by WMAP jllj : 



fi„ = 0.28 , (6) 

^DE = 0.72 , (7) 

Ho = 73 [km/s/Mpc] , (8) 

as = 0.74 . (9) 

For dark energy models, we set several equations of state. 

PDE = wpY,E, («; = -0.5,-0.8,-0.9,-1.0,-1.2), (10) 
where w = —1.0 corresponds to the cosmological constant. Even though w = —0.5 



has already been excluded and the strong energy condition is violated in the case of 
w = —1.2, we calculate these cases as well for comparison. 

According to the recent clustering results of XMM-Newton soft X-ray sources, the 
parameter ag depends on w The relation between ag and w can be flt by 

as = 0.34(±0.01)n-^(^'"''") , (11) 

7(fi„, w) = 0.22(±0.04) - 0A0{±0.05)w - 0.052(±0.040)n^ . (12) 

Using this relation, we compute as for several w cases (Tabled]). 

In order to avoid the divergence of the density fluctuation in the limit of large k 
and since the observed distribution is intrinsically discrete, it is necessary to consider 
the density fleld p{x] R) at the position x smoothed over some scale R, which is related 
to the unsmoothed density fleld p{x) as 

p{x;R)= I d^yWi\x-y\;R)p{y) 

/d'^h ~ 
^^W{kR)p{k)e-"'- , (13) 
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where W denotes a window function and W and p represent the Fourier transforms of the 
corresponding quantities. In this paper, after the calculations, we perform smoothing 
with a spherical top-hat window function, 

3(sinx — xcosx) , , 

W = ^ ^ , 14 

which yield the density fluctuation 6{x; R) at the position x smoothed over the scale 
R. For simplicity we use 6 to denote S{x; R) unless otherwise stated. 

It is also shown that for the case of smoothing with the spherical top-hat filter with 
radius R, in the weakly nonlinear region, the skewness and kurtosis of the density 
distribution can be obtained analytically for a given cosmology. For example, in 
Einstein-de Sitter Universe, they are given as 

60712 ^ 62d\og[a\R)] ^ 7 fd\og[a^{R)] 



1323 3 dlogR 3 V dlogR 

2 dHog[a\R)] , 
+ 3 dilogRy 

In our previous paper |l3j , we choose the smoothing scale i? as a rather large value 
{8h~^ Mpc), because Lagrangian linear perturbation cannot approximate the strongly 
nonlinear regime well because of shell-crossing or caustics formation. Since our present 
analysis is based on N-body simulations, we can reproduce small scale structure as well, 
which permits us to obtain meaningful results by choosing smaller values for R. 

Based on 6 calculated by the procedures above, we compute the statistical 
quantities, i.e., the dispersion, skewness and kurtosis for 6 time slices [z = 5, 4, 3, 2, 1, 0). 
In order to obtain the statistical quantities, we generate 10 samples for the primordial 
density fluctuations by COSMICS for each dark energy model. 

It is worth noting that we normalized the amplitude of the density spectrum at 
recombination era here, while the linear density spectrum is normalized by density 
dispersion at present era in other analyses like [221 [231 IMl [25] . Therefore in the nonlinear 
regime, our results seem different from theirs. However, as we show later, when w 
approaches to 0, the growth of the fluctuation is suppressed, which is consistent to 
them. 



5. Numerical results 

For comparison to the previous results, we set R = S/i^^Mpc in the comoving Eulerian 
coordinates at the present time (a = 1). In the previous paper, since we applied 
Lagrangian linear perturbation, we must examine the degree of the validity of the 
perturbative approach. 

(a) Validity of Lagrangian linear perturbation 
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Table 2. The difference of the dispersion, skewness, and kurtosis between Lagrangian 
linear perturbation and N-body simulations a± z = 2 {R — 8h~^ Mpc). 



w 


dispersion 


skewness 


kurtosis 


-0.5 


6.7% 


-23.4% 


-45.7% 


-0.8 


-1.6% 


-24.7% 


-47.8% 


-0.9 


-2.5% 


-24.8% 


-48.0% 


-1.0 


-3.0% 


-25.0% 


-48.2% 


-1.2 


-3.6% 


-25.1% 


-48.4% 



At first, we compare the dispersion of tlie density distribution obtained by 
Lagrangian linear perturbation and N-body simulations. Figure [1] shows the time 
evolution of the dispersions of the density distribution for w = —1.0. In the early 
stages, since the evolution is well described by the linear theory, the dispersions based 
on the two methods are almost the same. For example, the difference is less than 0.7% 
at z = 5. In our previous paper, we compare non-Gaussianity of the density distribution 
at z = 2 where the difference becomes a little larger, about 3%. 

Even though the difference between Lagrangian linear perturbation and N-body 
simulations is still small until z = 2 for the dispersion, it is not shown to hold also for 
non-Gaussianity of the density distribution. Figs. [2] and [3] show the time evolution 
of the skewness and kurtosis of the density distribution for w = —1, respectively 
and the difference is about 25% for the skewness and about 50% for the kurtosis. 
Therefore, it can be deduced that non-Gaussianity of the density distribution tend 
to be underestimated by Lagrangian linear perturbation, even for the region where the 
evolution of the dispersion is well described by it. 

Regardless of the fact above, the conclusion of our previous paper that the difference 
of non-Gaussianity between w = —1 and w = —0.8 or —1.2 is about 1% at z = 2 
remains unchanged, since the difference between Lagrangian linear perturbation and 
N-body simulations is also about 25% for the skewness and about 50% for the kurtosis 
for w = —0.8 or —1.2. We summarize the difference of the statistical quantities between 
Lagrangian linear perturbation and N-body simulations at z = 2 for various w (Table [2]). 

(b) Statistical quantities (R = 8h^^ Mpc) 

Next we compare the dispersion in several dark energy models based on N-body 
simulations. Now that we are using N-body simulations, we can go into the strongly 
nonlinear regime. Figure H] shows the time evolution of the dispersion of the density 
distribution in several dark energy models. The growth of the dispersion monotonously 
continues until z = 0. The tendency of the dispersion among dark energy models is 
similar to that obtained by Lagrangian linear perturbation. As the value of w becomes 
larger (approaches zero), the dispersion becomes smaller. This can be explained as 
follows: for the model with larger w, dark energy dominates at an earlier era, and 
the expansion of the Universe starts accelerating earlier, while it can be shown that 
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Figure 1. The dispersion of the density distribution based on Lagrangian linear 
perturbation and N-body simulations {R = 8/i~^Mpc) for w = —1. Lagrangian linear 
perturbation keeps to be a good approximation until z ^ 2. 

6 



5 

4 

I 3 
2 
1 


14 12 10 8 6 4 2 

z 

Figure 2. The skewness of the density distribution based on Lagrangian linear 
perturbation and N-body simulations [R — 8ft.""'^Mpc) for the case of w = —1. The 
difference of the skewness which is already about 5% at z = 15 increases. 

the growth of the density fluctuation is slower in the accelerating stage than in the 
matter-dominant stage. 

Figs. [5] and [6] show the time evolution of the skewness and kurtosis of the density 
distribution, respectively. It is different from that of the dispersion, in that the growth 
of them are not monotonous, especially after z = 1. Even though until z = 1, like 
the dispersion, these quantities are smaller as the value of w becomes larger, this order 
does not hold at the late time. Except the case w = —0.5 whose dispersion is not large 
enough, they are larger as the value of w becomes larger. 
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Figure 3. The Kurtosis of the density distribution based on Lagrangian linear 
perturbation and N-body simulations {R = Sh'^Mpc) for the case w — —1. The 
difference of the kurtosis which is already about 7% at z = 15 increases. 
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Figure 4. The dispersion of the density distribution in several dark energy models 
{R = 8h~-^Mpc). The growth of the dispersion continues monotonously. As the value 
of w becomes larger (approaches zero), the dispersion becomes smaller. The difference 
of the dispersion between w = — 1 and w = —0.8 is about 2.4% at z — 2, and become 
as large as 7.5% at z = 0. In the same way, the difference of the dispersion between 
w — —1 and w = —1.2 is 1.0% at z = 2 and 5.3% at z = 0. The difference between 
w = —1 and w = —0.5 is already 13% at z = 2. 



The detailed difference of the dispersion, skewness, and kurtosis for R = 8/i~^Mpc 
between w = —1.0 and other w at z = 2 and z = are summarized in tables [3] and HI 
respectively. 

Regardless of the change of the order of the w dependence, the differences of both 
the skewness and the kurtosis among w = —1.2, w = —1 and w = —0.8 are much less 
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Figure 5. The skewness of the density distribution in several dark energy models 
(i? = 8/i~^Mpc). The skewness oscillate between z = 1 and z = 0. We can easily 
distinguish the case for w = —0.5 from the results. The difference of the skewness 
between w = —1 and w — —0.8 is about 0.4% at z = 2, then spreads out about 
0.7% at 2; = 0. In the same way, the difference of the skewness between w = —1 and 
w = —1.2 is about 0.2% at z = 2, then spreads out about 0.3% at z = 0. 
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Figure 6. The kurtosis of the density distribution in several dark energy models 
{R = 8/i~^Mpc). The kurtosis oscillate between z = 1 and z = Q. As well as the 

skewness, the difference of the kurtosis between w = —1 and w = —0.8 is about 0.8% 
at 2 = 2, then spreads out about 1.2% at z = 0. In the same way, the difference of 
the kurtosis between w = —1 and w = —1.2 is about 0.3% at z = 2, then contract to 
about 0.2% at z = l. 

than !% ai z — 2 but turn out to become as large as 1% (skewness) and 2% (kurtosis) 
at ^ = 0. For the case o!w = —0.5, the values of both the skewness and the kurtosis are 
obviously less than that for the case of ty = —1.0 even ai z — 2 (about 2% and about 
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Table 3. The difference of the dispersion, skewness, and kurtosis between w — —1.0 
and other w at z = 2 (i? = S/i^^ Mpc). 



w 


dispersion 


skewness 


kurtosis 


-0.5 


-13.1% 


-2.1% 


-4.6% 


-0.8 


-2.4% 


-0.37% 


-0.80% 


-0.9 


-0.93% 


-0.14% 


-0.32% 


-1.2 


1.0% 


0.15% 


0.34% 



Table 4. The difference of the dispersion, skewness, and kurtosis between w = —1.0 
and other w a± z — {) {R — 8h^^ Mpc). 

w dispersion skewness kurtosis 

-0.5 -26.2% -0.09% 6i)% 

-0.8 -7.5% -0.67% -1.2% 

-0.9 -3.4% -0.24% -0.31% 

-1.2 5.3% 0.25% -0.19% 



5%, respectively). 

(c) Statistical quantities (R = 2h^^ Mpc) 

In this paper, since we apply N-body simulations, we can analyze smaller structure. 
In other words, even if we choose smaller value for R, we can obtain reliable results. 
Therefore, next we set the smoothing length R = 2h^^Mpc in the comoving Eulerian 
coordinates at the present time (a = 1). Even though the nonlinearity in this scale 
at z = is so strong that the corresponding galaxy distribution is influenced by the 
processes other than gravity such as the bias related with galaxy formation, it is worth 
investigating whether there exist other kind of scaling between higher order cumulants 
and the dispersion of the density distribution in this region. 

Figure [7] shows the time evolution of the dispersion of the density distribution in 
several dark energy models. The tendency for the growth of the dispersion is similar 
to that for R = 8/i^^Mpc. The dispersion monotonously grows until z = and if 
the equation of state of dark energy w increases, the dispersion would decrease. Since 
we choose smaller smoothing scale, the dispersion is larger than that in the case of 
R = 8/i~^Mpc for a fixed z and w. 

Figs. [H] and M show the time evolution of the skewness and kurtosis of the density 
distribution, respectively. As in the case of i? = 8/i~^Mpc, the growth of them is 
not monotonous. It is shown that after z = 1 for a fixed z, the skewness and 
kurtosis increases if the equation of state of dark energy w increases including the case 
w = —0.5. This suggests that in sufficiently nonlinear region, say a > 2, there exists 
some scaling between the higher order cumulants and the dispersion < 5" >cOC a^'-"'', 
where p{n) < 2n — 2, at least inside the scope of our calculations. 

The detailed difference of the dispersion, skewness, and kurtosis for R = 2h~^Mpc 
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Figure 7. The dispersion of the density distribution in several dark energy models 
[R = 2/i~-'^Mpc). The tendency of the difference of the dispersion is similar to that 
in the case of i? = 8/i~^Mpc. The difference of the dispersion between w = —1 and 
w — —0.8 is about 8.7% at z = 0. Similarly the difference of the dispersion between 
w — —1 and w = —1.2 is about 6.0% at z = 0. The difference of the dispersion between 
w = —1 and w = —0.5 is quite large (about 17% already at z = 2). 



Table 5. The difference of the dispersion, skewness, and kurtosis between w = —1.0 
and other w ai z^2 {R^ 2h^^ Mpc). 



w 


dispersion 


skewness 


kurtosis 


-0.5 


-16.8% 


-7.8% 


-13.0% 


-0.8 


-3.1% 


-1.2% 


-1.9% 


-0.9 


-1.2% 


-0.48% 


-0.81% 


-1.2 


1.3% 


0.53% 


0.97% 



between w = —1.0 and other w at z = 2 and z = are summarized in tables [5] and El 
respectively. 

Regardless of the change of the order of the w dependence, the differences of both 
the skewness and the kurtosis among w = —1.2, w = —1 and w = —0.8 are much less 
than 2% at z = 2 but turn out to become as large as 3% (skewness) and 9% (kurtosis) 
at z = 0. For the case of w = —0.5, the values of both the skewness and the kurtosis are 
obviously less than that for the case of w = —1.0 even at 2; = 2 (about 8% and about 
13%, respectively). 

(d) Probability distribution function 

In order to clarify the difference between the case of w = —0.8 and w = —1.0, 
we present raw data of PDF of the density distribution at 2; = for R = 8/i^^Mpc 
(Figure [TO!) . and R = 2/i~^Mpc (Figure [TTj) . respectively. When the nonlinear effect 
promotes the formation of the dense structure and restrains the evolution of the voids. 
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Figure 8. The skewness of the density distribution in several dark energy models 
(i? = 2/i~^Mpc). The growth of the skewness stops after z = 1. The difference of 
the skewness between w = —1 and w — —0.8 is about 3.0% at z = 0. At some time 
between z — 1 and z = 0, the skewness of w = —0.8 becomes larger than that of 
w = —1.0. In the same way, the difference between w — ~l and w = —1.2 is about 
2.4% at z = 0. At some time between z = 1 and z = 0, the skewness oi w — —1.2 
becomes smaller than that of it; —1.0. The dilference between w = —1 and w = —0.5 
is as large as 3.4% even at z = 0. 

Table 6. The difference of the dispersion, skewness, and kurtosis between w — —1.0 
and other w at z = (i? = 2h^^ Mpc). 



the degree of asymmetry of PDF (skewness) grows obviously. On the other hand, if 
the edge of PDF spreads widely without reduction of peak, the "peakiness" of PDF 
(kurtosis) grows evidently. 

(e) Error analysis 

Until now, we provided our results without error bars, even though they are 
obtained by averaging 10 samples for each dark energy model. However, in order to 
check whether these results can be used directly as a cosmological tool, the statistical 
error analysis is necessary. Here we estimate the variances among samples for the case 
of w = —1.0 by increasing the number of samples up to 50. The results are showed in 
Figs. [T2lfT^ As a result, contrary to our expectation, the variance is much larger than the 
difference between the difference of the non-Gaussianity between w = —1 and the other 
cases. Furthermore, it turns out that although the variance of the dispersion can be 



w 



dispersion skewness kurtosis 



-0.5 
-0.8 
-0.9 
-1.2 



-32.2% 3.4% 7.7% 

-8.7% 3.0% 4.5% 

-3.9% 1.2% 2.6% 

6.0% -2.4% -9.0% 
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Figure 9. The kurtosis of the density distribution in several dark energy models 
[R = 2/i~^Mpc). The growth of the kurtosis stops at z = 1. The difference of the 
kurtosis between w = —\ and w = —0.8 is about 4.5% at z = 0. At some time between 
z = 1 and z = 0, the kurtosis of w = —0.8 becomes larger than that oi w = —1.0. In 
the same way, the difference between w = —1 and w = —1.2 is about 9% at z = 0. 
At some time around z — 2, the kurtosis of w = —1.2 becomes smaller than that of 
w = —1.0. The difference between w = —1 and w — —0.5 is as large as 8% even at 
z = 0. 
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Figure 10. The PDF of the density distribution in w ~ —0.8 and w = —1.0 at z = 
{R = 8h^^ Mpc). Because of smoothing for the PDF, we pick up the PDF with low 
resolution {AS = 2). Since the PDF of w = —1.0 spreads more widely than that of 
w — —0.8 for almost the same dispersion, the non-Gaussianity of the PDF of ui = —1.0 
is larger than that of w = —0.8. 

decreased , the variance of the non-Gaussianity is not necessarily improved by increasing 
the number of samples. The detailed results are shown in Figs. [T2lfT^ Therefore, for the 
practical use of this analysis to constrain the dark energy models, further improvements 
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Figure 11. The PDF of the density distribution in w — —0.8 and w = —1.0 at 2; = 
{R = 2h~^ Mpc). Because of smoothing for the PDF, we pick up the PDF with low 

resolution (A5 = 4). Even though the PDF of u> = —1.0 spreads more widely than 
that of w = —0.8, since the dispersion of the PDF of w = —1.0 is larger than that of 
w = —0.8, the resulting non-Gaussianity becomes smaller for the case oi w = —1.0. 
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Figure 12. The dispersion of the density distribution in the case ot w = —1.0 
{R = 2/i~^Mpc). Here we analyze the variance of samples. When we increase number 
of samples, the variance slightly decreases. The variance is about 7% at 2; = in both 
cases. 



for the calculation of the statistical quantities, which seems to be worth considering is 
necessary. 
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Figure 13. The skewness of the density distribution in the case oi w = —1.0 
[R — Here we analyze the variance of samples. Even if we increase number 

of samples, the variance of the non-Gaussianity is not always improved. The variance 
is about 18% at z = 0. 
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Figure 14. The kurtosis of the density distribution in the case oi w = —1.0 
(i? = 2/i~^Mpc). Here we analyze the variance of samples. Even if we increase number 
of samples, the variance of the non-Gaussianity is not always improved. The variance 
is about 46% at z = 0. 

6. Summary 

To clarify the nature of dark energy is, without doubt, one of the most important tasks in 
modern cosmology. From the phenomenological viewpoint, it is important to constrain 
the effective equation of state of dark energy w by observations. Even though we have 
obtained the constraint as w < —0.90 (95 % confidence limit assuming w > — 1) [llj by 
combining WMAP data with other astronomical data, to pin down the value of ly, it is 
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necessary to propose other independent observational tools. 

For this purpose, we consider non-Gaussianity of the density distribution which 
is generated by the non-linear dynamics. In this paper, as a natural extension to our 
previous work [13] where we follow a semi-analytic approach based on Lagrangian linear 
perturbation theory, we apply N-body simulations based on P^M codes. In terms of 
the density fluctuation obtained by N-body simulations, we compute the skewness and 
kurtosis of the smoothed density distribution with the appropriate definition in that they 
become constant in extremely weakly nonlinear stage. They are important statistical 
quantities denoting the deviation from the Gaussian distribution. By considering several 
constant-ty dark energy models {w = —0.5, —0.8, —0.9, —1.0, —1.2), we present the in- 
dependence of the time evolution of non-Gaussianity of the density distribution. 

First we compare the results for the weakly nonlinear region at 2; = 2 obtained by N- 
body simulations and Lagrangian linear perturbation by choosing the same smoothing 
length R = Sh'^ Mpc. In this region where non-Gaussianity become smaller as w 
becomes larger, we regarded the results by Lagrangian linear perturbation as reliable 
in the previous work. As a result, even though the dispersion obtained by Lagrangian 
linear perturbation agrees well with the ones obtained by N-body simulations until 
2; ~ 2, non-Gaussianity tend to be underestimated by Lagrangian linear perturbation. 

This result is attributed to the fact that second order Lagrangian perturbation 
is necessary to express second order Eulerian perturbation [371 EH], mathematically 
and the peculiar acceleration perpendicular to the peculiar velocity can not be picked 
up by Lagrangian linear perturbation, physically. Therefore, we learn the fact that 
"Lagrangian linear perturbation is not necessarily better than second order Eulerian 
perturbation". However, regardless of the underestimation of non-Gaussianity, the 
conclusions of the previous paper that the relative difference of non-Gaussianity between 
w = —1.0 and w = —0.8 ot w = —1.2 is about 1% at z = 2 for R = 8h~'^Mpc remains 
unchanged because the degree of the underestimation in Lagrangian linear perturbation 
is almost w-independent. 

Next, we analyze the statistical quantities in the strongly nonlinear region where the 
previous perturbative approach is obviously invalid. We find that the relative difference 
of non-Gaussianity among several constant-w models spread as 2 — > 0, even though 
the growth of non-Gaussianity is not monotonous in that region and the order changes 
around z = 0.5 except w = —0.5 model whose dispersion is not so large enough. For 
example, at 2 = the relative difference between non-Gaussianity of w = —0.8 model 
and that of w = —1.0 model is 0.67% for the skewness and 1.2% for the kurtosis, while 
at 2; = 2 they are —0.52% for the skewness and —1.3% for the kurtosis. Similarly, at 
z = the difference between non-Gaussianity of w = —1.2 model and that of w = —1.0 
model is 0.25% for the skewness and —0.19% for the kurtosis, while at 2; = 2 they are 
0.15% for the skewness and 0.34% for the kurtosis. 

Making use of N-body simulations, we also analyze the smaller structure by setting 
the smoothing scale R = 2/i~^Mpc, even though the nonlinearity in this scale at 2 = 
is so strong that the galaxy distribution is influenced by the process other than gravity. 
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For this case, since we look at the smaller scale, the dispersion and non-Gaussianity 
of the density distribution is larger than that in the case of R = 8h~^Mpc for a fixed 
z and w. Relative difference of non-Gaussianity among several constant-w models are 
also greater than the case of i? = 8/i~^Mpc. 

For example, at z = the relative difference between non-Gaussianity ofw = —0.8 
model and that of w = —1.0 model is 3.0% for the skewness and 4.5% for the kurtosis. 
Similarly, the relative difference between the non-Gaussianity of w = —1.2 model 
and that of w = —1.0 model is —2.4% for the skewness and —9.0% for the kurtosis. 
Furthermore, on the contrary to the weakly nonlinear region, for a fixed z we confirm 
non-Gaussianity become larger as w becomes larger after z = 0.5, including the case of 
w = -0.5. 

This suggests that in the sufficiently nonlinear region, say a > 2, there exists some 
scaling between the higher order cumulants and the dispersion < 5" >cOC cr^*^"\ where 
p{n) < 2n — 2, at least inside the scope of our calculations. Academically, not limited to 
the constraint for the dark energy models, this result seems to worth investigating related 
with the nature of the nonlinear gravitational clustering in the expanding Universe. 

It is necessary to mention the possibility of examining our results. Recently, several 
galaxy redshift survey projects have been progressing [391 SOI SI]- In these projects, 
many galaxies within a region {z < 0.3 for 2dF, z < 0.5 for SDSS, and 0.7 < z < 1.4 for 
DEEP2) have been observed. These projects show the latest distribution of galaxies, 
which form strongly nonlinear structure. Actually, the skewness and kurtosis of the 
distribution of SDSS Early Data Release had been computed [12]. The results suggested 
the SDSS imaging data can enable to determine the skewness and kurtosis up to 1% and 
less than 10%, respectively. Therefore, in future projects we can expect to determine 
the skewness and kurtosis with high accuracy. However, based on our error analysis, 
the variance of the non-Gaussianity turned out to be about 50% for 10 samples and it 
does not increased even by increasing the number of samples up to 50. This suggests 
that our analysis, by itself, is still immature and further improvements are necessary for 
the practical use to constrain the property of dark energy. 

As mentioned above, for smaller scale like R = 2h~^Mpc now, the complex 
physical processes about galaxy formation affect the galaxy distribution. In order to 
use observational data such as the PSCz survey [43] which is an all sky survey at the 
present epoch, it is necessary to create mock PSCz catalogs from our results making 
use of the PSCz selection function [44J. Even though this is from excessively academical 
interest, if-dependence of the non-Gaussianity of the actual galaxy distribution seems 
to be worth investigating. On the other hand, it is known that weak lensing surveys 
can potentially provide us with precision maps of the projected density up to redshifts 
around 1 [l5l HBl HTJ HH] • Even though we need another step of obtaining the convergence 
field which can be written as the projection of the matter density along the line of sight, 
the skewness and kurtosis of the convergence field can be tested by weak lensing surveys. 
Obviously, our results obtained in this paper help to derive these quantities. 

Finally, it is worth noting that primordial non-Gaussianity which is generated by 
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inflation [2P|. From WMAP observations, a limit for the non-linear coupling parameter 
has been established |19] which does not deny the existence of non-Gaussianity in the 
primordial density fluctuations. Even though we cannot disentangle them generally, 
we can calculate non-Gaussianity of large scale structure for a given non-Gaussian 
initial condition and compare the results with the one obtained by the Gaussian initial 
condition. This requires further investigations and we hope to report results in a separate 
publication. 
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